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SUMMARY 


The usefulness of NASTRAN for nuclear applications has been investigated. 
For this purpose the extent to which suitable solutions may be obtained for one 
physics problem and two engineering type problems is traced. NASTRAN appears 
to be a practical tool to solve one-group steady-state neutron diffusion equa- 
tions. Transient diffusion analysis may be performed after new levels that 
allow time-dependent temperature calculations are developed. 

NASTRAN Piecewise Linear Analysis may be applied to solve those plasticity 
problems for which a smooth stress-strain curve can be used to describe the non- 
linear material behavior. The accuracy decreases when sharp transitions in the 
stress-strain relations are involved. Improved NASTRAN usefulness will be 
obtained when nonlinear material capabilities are extended to axisymmetric ele- 
ments and to include provisions for time-dependent material properties and creep 
analysis. Rigid Formats 3 and 5 proved to be very convenient for the buckling 
and normal -mode analysis of a nuclear fuel element. 


INTRODUCTION 


In addition to NASTRAN finding its way into leading U.S. aerospace indus- 
tries and research institutes, it is also more and more used in other areas of 
structural analysis such as civil engineering and the automotive industry. In 
the past much duplication has occurred with respect to the development of finite- 
element computer programs, but at present it is worthwhile to investigate the 
applicability of NASTRAN before deciding to write a finite -element routine or to 
add large extensions to existing programs. NASTRAN is an easy-to-use general- 
purpose program with assured continuous management and maintenance and is avail- 
able at a low price. 

The domain of nuclear problems is an important field of application outside 
the area for which NASTRAN was originally developed, and the present paper dis- 
cusses the use of NASTRAN for some characteristic nuclear engineering and physi- 
cal problems. The recently released Level 15 allows the solution of physical 
problems since a functional for neutron diffusion calculations, which shows a 
close relationship to one used for thermal analysis, was derived. Except for 
the computation of a neutron-flux distribution, the paper deals with two struc- 
tural example problems j one discusses the elastic and plastic deformation of 
nuclear fuel cladding, and the other concerns the normal-mode and buckling ana- 
lysis of a fuel element. 
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The example problems chosen do not cover the whole area of finite- element 
applications suitable for nuclear research and design nor are all details of 
them investigated. Therefore, the paper cannot provide sufficient information 
needed to decide whether to use NAS TRAN or not and may only serve as a guideline 
to more directed research. 


ONE-GROUP NEUTRON DIFFUSION 


Recently, different authors (refs. 1, 2, and 3) have applied finite -element 
techniques to solve problems from nuclear physics. The equations used to 
describe heat transfer and neutron diffusion show a close similarity, and only 
a few minor modifications of heat transfer programs are necessary to allow 
neutron diffusion calculations. Functionals for both phenomena will be discussed 
and NASTRAN alterations to give neutron-flux distributions will be outlined. 

The results for one particular application are presented. 

For the case in which all neutrons are assumed to belong to the same energy 
(velocity) group, the equation of neutron diffusion reads (ref. 4) 


_ i M + d — c?0 + S = 0 

v at a Xi a Xi 


(1> 


Here the summation convention is used while v stands for velocity, D for the 
diffusion coefficient, 0 for the scalar flux, a for the macroscopic absorp- 
tion cross section, and source term S for the rate of production of neutrons 
per unit volume per second. 

The boundary conditions applied to equation (l) involve continuity of flux 
and current at interfaces, and vanishing of flux at extrapolated boundaries or 
return current at the actual boundaries. The general form of the boundary con- 
ditions reads 


P + 70 + q = 0 (2) 

bn 


where y and q are constants. Equations (l) and (2) can be combined to obtain 
the functional: 
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where the triple integrals are to be extended over the volume V, and the double 
integrals over the surface A of this volume. Functional F reaches a m-iniTmim 
value with respect to all admissible variations of flux distribution 0. 
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The functional most commonly used for finite -element heat conduction compu- 
tation is 


P = 



- ae ro e 


| 0 2 ^dA 


00 


which reaches a minimum value with respect to all admissible variations of tem- 
perature distribution 9. The symbol Q° is used for volumetric heat genera- 
tion rate, p for density, and c for specific heat while A stands for the 
coefficient of heat conductivity. The coefficient of heat transfer between the 
structure and its surrounding medium is represented by a and the temperature 
of the surrounding medium by 0 m . Comparison of expressions (3) and (4) shows 
that duality of finite-element temperature programs is obtained when for neutron 


diffusion calculations the following substitutions are made: — for p j 

n V G 

D for A; S for Q°; - — for 8 m ; and y for a. Moreover, the heat 


capacity matrix has to be generated using a coefficient of thermal capacity 
equal to a and added to the heat conduction matrix. As yet NAS TRAN heat flow 
capability can be used for stationary problems only and in this way applications 
to the calculation of neutron-flux distributions have to be restricted to steady- 
state conditions. The calculations can be performed by carrying out the substi- 
tutions mentioned and the insertion of an ALTER package (fig. 1). In Level 15- 1.1 
the heat capacity matrix cannot be generated directly. However MAT1, MAT2, or 
MAT3 cards can be used to submit values for the mass density equal to those 
desired for the thermal capacity and module SMA2 can be executed to obtain the 
mass matrix. The ALTER package allows the generation of matrix MGG and the addi- 
tion of it to the heat conduction matrix KGGX, Module SCE2 removes additional 
matrix elements implied by MGG and undesired in the heat capacity matrix. The 
method described has been used to determine the flux distribution in the basic 
area of symmetry in a fuel-moderator mixture located between a square array of 
cruciform control rods. The rods are black to thermal neutrons and a uniform 
slowing density has been assumed. Figure 2 shows the geometry together with the 
idealization and calculated lines of constant neutron flux. The problem has also 
been investigated by Semenza et al, (ref. 2) and the agreement looks reasonable. 


DEFORMATION OF NUCLEAR-FUEL CLADDING 


One of the most important problems in the design of nuclear fuel elements 
is the interactive deformation of fuel pellets and the surrounding cladding 
(fig. 3). For most of the current cladding and fuel materials, the temperatures 
at which the loadings are applied will cause creep and plasticity to occur. 
Moreover, the nuclear irradiation leads to swelling, introduced by helium bubble 
formation, and embrittlement of the materials ; thus, two more complications must 
be dealt with in analyzing the deformations. 

With the exception of contact forces implied by the cladding, the pellet is 
only loaded thermally. The cladding is subject to 
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(1) Differential pressure between coolant on the outside and released 

fission gas on the inside 

(2) Contact pressure and axial forces due to swelling and thermal defor- 

mation of the fuel 

(3) Thermal load caused by heat transfer through the clad wall 

The magnitude of these loadings depends on time, neutron dose, and power output 
of the reactor. 

NASTRAN permits the calculation of deformations both elastically and 
plastically. A static solution could be obtained by use of solid-of- revolution 
elements. After a choice has been made with respect to the grid points where 
contact occurs, multipoint- constraint equations could provide expressions for 
the displacements of internal cladding grid points in terms of external pellet 
grid points and the undeformed gap. The idealization and the deformed shape are 
given in figure 4, and the results show close agreement with solutions obtained 
from other programs (ref. 5). 

Plasticity calculations of cladding deformations as a result of mechanical 
loading have been performed using Piecewise Linear Analysis. Because this 
feature is as yet not applicable for solid- of -revolution elements, the cladding 
wall has been idealized by parallel plate elements (fig. 5)- I* 1 using this 

model, all axial displacements were restrained and all elements were kept par- 
allel during deformation while interelement distances remained constant. More- 
over all collinear grid points remained collinear and on the same radii. All 
these restrictions could easily be accomplished by multipoint and single point 
constraints leading to a number of independent degrees of freedom equal to 1. 

The application of this model means that radial stresses are neglected. These, 
however, proved to be of minor importance and are partly counterbalanced by the 
assumption of constant interelement distances. 


During each pass through the Piecewise Linear Analysis loop a new value for 
the stress-strain slope is estimated by using 


E 1 = 



(5) 


where and e * are the uniaxial equivalences of stress and strain obtained 

after application of the ith load increment and 




+ 


*°i+i 

Aa^ 




( 6 ) 


i+1 = 
a 



( 7 ) 
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where Aa^ denotes the load increment during the ith load factor and F denotes 
a user-supplied tabular function of stresses and strains. 

If a strain-hardening law is applied, the model requires that within a 
small increase of load, all elements deform beyond the elastic limit. The over- 
all rigidity decreases considerably in the plastic stage and a load increment 
beyond the yielding point causes a strain increment which may be many times 
greater than that registered during the elastic deformation. 

This is why the use of expression (6) to predict the strain due to the next 
load factor can lead to greatly underestimating the strain and, in such cases, 
the value used for E 1 will be much too large. The present levels of NASTRAN 
do not contain a check on the accuracy of E^. For the model described, accept- 
able results could only be obtained by choosing the piecewise linear factors in 
such a way that each of the nonlinear elements will be loaded up to the yielding 
point for one particular load factor. 

This approach is of course very laborious and not feasible for structures 
with many nonlinear elements. However, the accuracy can be improved if the 
relationship between stresses and strains can be described by a smooth curve 
and if small load increments are used. 


Thermal stresses are an important factor in the analysis of the cladding. 
When these stresses alternate, as is the case during startups and shutdowns of 
the reactor, they become of even more importance and a stepwise increase of the 
cladding diameter, called "thermal ratchetting, " may be observed. 

In Piecewise Linear Analysis (PLA) , thermal loading can be taken into 
account by generating the relevant load vector separately and adding it stepwise 
to the mechanical loading after the mechanical loading has been fully accumu- 
lated during repeated execution of the PLA loop. The deviating values of E* 
and the alternate direction of strain increments for elements compressed by the 
thermal strain could be accomplished by application of Direct Table Input. 

These tables can substitute the regular element summary, connection, and prop- 
erty tables for nonlinear elements, but once again the approach proved to be 
laborious and useless for practical applications. NASTRAN did not allow creep 
calculations or the consideration of time -dependent material properties intro- 
duced by the influence of nuclear irradiation. 


BUCKLING AND NOEMAL-MODE ANALYSIS OF FUEL ELEMENT 


In many types of nuclear reactors, the heat generated by nuclear fission 
is removed by a gaseous or liquid coolant flowing longitudinally along the fuel 
rods. Considering the buckling of these rods when subjected to the flowing 
environment, it does not matter whether the fluid is assumed to flow through or 
along the tube, and the critical value of the fluid velocity can be found by 
solving the differential equation for steady flow of the fluid through the tube: 
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El ■* pv 2 ti. = p 
cix^ dy? 

Here z stands for the transverse displacement, v for fluid velocity, and El 
for the bending stiffness of the tube. This agrees (ref. 6) with the equation 
for a pipe carrying a lateral load p and an axial load of pv2. When Z is 
used for the length it means that buckling of the pipe will occur when 

<"* . El(2) 2 

therefore, the critical velocity is 


v 


3t 

Z 



NASTRAN allows accurate 1 determination of the buckling load, and therefore the 
critical velocity, when the pipe is idealized by only a few collinear BAR ele- 
ments. Therefore, this program may also be used for the computation of crit- 
ical velocities in complete fuel elements. In order to reduce the number of 
degrees of freedom, a NASTRAN fuel-element analysis has been performed by use of 
an artificial element with fewer fuel rods than usual but with a length-to-width 
ratio in agreement with existing fast-breeder reactor designs. The element con- 
sidered contains 37 fuel rods arranged in a hexagonal pattern with an 8-mm pitch. 
The 6-mm-diameter rods occupy the full 800-mm length of the element. 

To support the rods in the transverse directions, the element is equipped 
with four intermediate honeycomb-like grids; thus, all transverse displacements 
are equated at grid locations. The element is assumed to be clamped in the 
lower element supporting plate, whereas only longitudinal displacements are pos- 
sible in the upper plate. 

Every fuel rod is idealized by 10 collinear BAR elements of equal length 
while each side of the hexagonal prismatic container is modeled by a rectangular 
pattern of QUAD1 elements with 10 elements over the length and four elements 
over the width. Figure 6 gives the undeformed element together with normal 
modes 1, 2, and 5 . The first four modes show vibrations of the fuel rods and 
mode 5 represents the impact of the lowest container frequency on the behavior 
of the fuel element also. When liquid sodium is used for cooling, the critical 
velocity is found to be 5*7 m/sec. Matrix partitioning has been applied to 
obtain 195 degrees of freedom in the a- set. For the problem concerned, this 
figure is too high to produce a significant reduction in CRJ time. Checkpointing 
with restarting proved to be profitable when other eigenvalue ranges of interest 
had to be investigated. 


1 


0.16 percent deviation for three elements. 
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CONCLUDING REMARKS 


The use of NASTRAN for the solution of three nuclear problems has been 
investigated. In view of the fact that NASTRAN was developed for solution of 
problems in aerospace structures, it is not astonishing that it Is somewhat 
limited in the scope of nuclear problems that it will solve and requires more 
user effort than special-purpose nuclear programs. In summary, the conclusions 
from the analysis of the chosen examples can be stated as follows: 

1. Steady-state one-group neutron-flux distributions can be computed very 
conveniently by the use of NASTRAN. One of the next official levels featuring 
transient-temperature calculations will also permit the solution of time- 
dependent diffusion equations. The applicability to multigroup diffusion and 
other problems from nuclear physics may be investigated. 

2. At this time, Rigid Format 6 seems unsuitable when a strain-hardening 
law is used and alternating thermal loads are involved. In order to improve the 
present approach used in Piecewise Linear Analysis (PLA) it may be advisable to 
repeat the calculation in those cases where a large deviation between the esti- 
mated and actual strains occurs. Reference 7> or 9 may be useful if an alter- 
native approach is shown to be desirable. Other areas of interest for nuclear 
users will be to have PLA available for solid-of-revolution elements also and to 
be able to apply NASTRAN for creep calculations. 

3. Buckling and normal-mode analysis can be applied directly to nuclear 
problems and to permit calculations of critical fluid velocities. The use of 
matrix partitioning reduces the number of independent degrees of freedom without 
loss of accuracy. More emphasis on creep applications is also desirable for 
these rigid formats, and a provision to perform creep-buckling analysis would 
mean a considerable improvement. 

The criticism expressed did not show to full advantage the many important fea- 
tures such as checkpointing, substructuring, and multipoint constraints which 
are not found in any known computer code for nuclear use. These features give 
NASTRAN the lead over its competitors. The general and flexible applicability 
of NASTRAN constitutes its main power, and especially when more effort can be 
spent to improve nonlinear material capabilities, NASTRAN will have a large 
potential in the area of nuclear research and design. 
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Figure 1.- Executive Control Deck for neutron diffusion calculations. 



Figure 2.- Idealization and lines of constant neutron flow in a fuel-moderator 
mixture area between cruciform control rods. 
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Figure 6.- Normal modes of nuclear fuel element. 




